In this report, we reproduce the analyses using data from fMRI study 1 reported in Supplementary Material.

prep data

First, we load the relevant packages, define functions and plotting aesthetics, and load and tidy the data.

load packages

library(pacman)
pacman::p_load(tidyverse, purrr, fs, knitr, lmerTest, ggeffects, kableExtra, boot, devtools, install = TRUE)
devtools::install_github("hadley/emo")

define functions

source("https://gist.githubusercontent.com/benmarwick/2a1bb0133ff568cbe28d/raw/fb53bd97121f7f9ce947837ef1a4c65a73bffb3f/geom_flat_violin.R")

# MLM results table function
table_model = function(model_data) {
  model_data %>%
  broom.mixed::tidy(conf.int = TRUE) %>%
  filter(effect == "fixed") %>%
  rename("SE" = std.error,
         "t" = statistic,
         "p" = p.value) %>%
  select(-group, -effect) %>%
  mutate_at(vars(-contains("term"), -contains("p")), round, 2) %>%
  mutate(term = gsub("cond", "", term),
         term = gsub("\\(Intercept\\)", "intercept", term),
         term = gsub("condother", "other", term),
         term = gsub("condself", "self", term),
         term = gsub("topichealth", "topic (health)", term),
         term = gsub("self_referential", "self-referential", term),
         term = gsub("self_relevance", "self-relevance", term),
         term = gsub("social_relevance", "social relevance", term),
         term = gsub(":", " x ", term),
         p = ifelse(p < .001, "< .001",
             ifelse(p == 1, "1.000", gsub("0.(.*)", ".\\1", sprintf("%.3f", p)))),
         `b [95% CI]` = sprintf("%.2f [%0.2f, %.2f]", estimate, conf.low, conf.high)) %>%
  select(term, `b [95% CI]`, df, t, p)
}

simple_slopes = function(model, var, moderator, continuous = TRUE) {
  
  if (isTRUE(continuous)) {
    emmeans::emtrends(model, as.formula(paste("~", moderator)), var = var) %>%
      data.frame() %>%
      rename("trend" = 2) %>%
      mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", trend, asymp.LCL, asymp.UCL)) %>%
      select(!!moderator, `b [95% CI]`) %>%
      kable()  %>%
      kableExtra::kable_styling()
    
  } else {
    confint(emmeans::contrast(emmeans::emmeans(model, as.formula(paste("~", var, "|", moderator))), "revpairwise", by = moderator, adjust = "none")) %>%
      data.frame() %>%
      filter(grepl("control", contrast)) %>%
      mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", estimate, asymp.LCL, asymp.UCL)) %>%
      select(contrast, !!moderator, `b [95% CI]`) %>%
      arrange(contrast) %>%
      kable()  %>%
      kableExtra::kable_styling()
  }
}

define aesthetics

palette_condition = c("self" = "#ee9b00",
                      "control" = "#bb3e03",
                      "other" = "#005f73")

palette_sharing = c("#0a9396", "#ee9b00")
palette_roi = c("self-referential" = "#ee9b00",
               "mentalizing" = "#005f73")
palette_dv = c("self-relevance" = "#ee9b00",
               "social relevance" = "#005f73",
               "sharing" = "#56282D")
palette_topic = c("climate" = "#E6805E",
                 "health" = "#3A3357")

plot_aes = theme_minimal() +
  theme(legend.position = "top",
        legend.text = element_text(size = 12),
        text = element_text(size = 16, family = "Futura Medium"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        axis.text = element_text(color = "black"),
        axis.line = element_line(colour = "black"),
        axis.ticks.y = element_blank())

load and tidy data

merged_all = read.csv("../data/study1_data.csv")

merged = merged_all %>%
  filter(outlier == "no" | is.na(outlier)) %>%
  group_by(pID, atlas) %>%
  mutate(parameter_estimate_std = parameter_estimate / sd(parameter_estimate, na.rm = TRUE)) 

merged_wide = merged %>%
  filter(atlas %in% c("self-referential", "mentalizing")) %>%
  select(site, pID, trial, topic, cond, value, self_relevance, social_relevance, atlas, parameter_estimate_std) %>%
  spread(atlas, parameter_estimate_std) %>%
  rename("self_referential" = `self-referential`)

merged_wide_alt = merged %>%
  filter(atlas %in% c("pnas_self", "pnas_mentalizing_nopc")) %>%
  select(site, pID, trial, topic, cond, value, self_relevance, social_relevance, atlas, parameter_estimate_std) %>%
  spread(atlas, parameter_estimate_std) %>%
  rename("self_referential" = pnas_self,
         "mentalizing" = pnas_mentalizing_nopc) 



sensitivity analyses

Given the high correlation between the preregistered Neurosynth ROIs, we conducted sensitivity analyses using ROIs from Scholz et al. (2017) A neural model of valuation and information virality.

In order to maximize the differentiation between the self-referential and mentalizing ROIs, we removed the PCC/precuneus cluster from the mentalizing ROI as it overlapped with the self-referential ROI.

ROI correlations

Compared to the preregistered Neurosynth ROIs (r = .94, 95% CI [.93, .94]), the correlation between the alternative ROIs are substantially reduced.

merged_wide_alt %>%
  rmcorr::rmcorr(as.factor(pID), mentalizing, self_referential, data = .)
## 
## Repeated measures correlation
## 
## r
## 0.5593866
## 
## degrees of freedom
## 5930
## 
## p-value
## 0
## 
## 95% confidence interval
## 0.5416457 0.5766294

H1

Is greater activity in the ROIs associated with higher self and social relevance ratings?

self-referential ROI

✅ H1a: Greater activity in the self-referential ROI will be associated with higher self-relevance ratings

mod_h1a =  lmer(self_relevance ~ self_referential + (1 + self_referential | pID),
               data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h1a = table_model(mod_h1a)

table_h1a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.56 [2.48, 2.63] 84.66 65.05 < .001
self-referential 0.03 [0.01, 0.06] 83.65 2.44 .017

summary

summary(mod_h1a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: self_relevance ~ self_referential + (1 + self_referential | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16765.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.4710 -0.7073  0.1362  0.6758  2.4472 
## 
## Random effects:
##  Groups   Name             Variance Std.Dev. Corr 
##  pID      (Intercept)      0.116819 0.34179       
##           self_referential 0.001551 0.03938  -0.81
##  Residual                  0.918352 0.95831       
## Number of obs: 6016, groups:  pID, 85
## 
## Fixed effects:
##                  Estimate Std. Error       df t value            Pr(>|t|)    
## (Intercept)       2.55656    0.03930 84.66264  65.045 <0.0000000000000002 ***
## self_referential  0.03178    0.01300 83.65164   2.444              0.0166 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## self_rfrntl -0.346

mentalizing ROI

✅ H1b: Greater activity in the mentalizing ROI will be associated with higher social relevance ratings

mod_h1b = lmer(social_relevance ~ mentalizing + (1 + mentalizing | pID),
               data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h1b = table_model(mod_h1b)

table_h1b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.66 [2.58, 2.75] 84.14 64.14 < .001
mentalizing 0.04 [0.01, 0.06] 82.92 3.10 .003

summary

summary(mod_h1b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: social_relevance ~ mentalizing + (1 + mentalizing | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 15861.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7979 -0.7199  0.1741  0.6482  2.6749 
## 
## Random effects:
##  Groups   Name        Variance   Std.Dev. Corr 
##  pID      (Intercept) 0.13456626 0.366833      
##           mentalizing 0.00005386 0.007339 -0.21
##  Residual             0.78577252 0.886438      
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##             Estimate Std. Error       df t value             Pr(>|t|)    
## (Intercept)  2.66386    0.04153 84.13768   64.14 < 0.0000000000000002 ***
## mentalizing  0.03545    0.01143 82.92434    3.10              0.00264 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## mentalizing -0.093

combined plot

predicted = ggeffects::ggpredict(mod_h1a, c("self_referential [-4.5:5]")) %>%
  data.frame() %>%
  mutate(roi = "self-referential",
         variable = "self-relevance") %>%
  bind_rows(ggeffects::ggpredict(mod_h1b, c("mentalizing [-4.5:5]")) %>%
              data.frame() %>%
              mutate(roi = "mentalizing",
                     variable = "social relevance"))

ind_data = merged_wide %>%
  select(pID, trial, contains("relevance"), mentalizing, self_referential) %>%
  rename("self-referential" = self_referential) %>%
  gather(variable, predicted, contains("relevance")) %>%
  mutate(variable = gsub("self_relevance", "self-relevance", variable),
         variable = gsub("social_relevance", "social relevance", variable)) %>%
  gather(roi, x, mentalizing, `self-referential`) %>%
  filter(!(variable == "self-relevance" & roi == "mentalizing") & ! (variable == "social relevance" & roi == "self-referential"))

(plot_h1 = predicted %>%
  ggplot(aes(x, predicted)) +
  stat_smooth(data = ind_data, aes(group = pID, color = roi), geom ='line', method = "lm", alpha = .1, size = 1, se = FALSE) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high, fill = roi), alpha = .3, color = NA) +
  geom_line(aes(color = roi), size = 2) +
  facet_grid(~variable) +
  scale_color_manual(name = "", values = palette_roi, guide = FALSE) +
  scale_fill_manual(name = "", values = palette_roi, guide = FALSE) +
  labs(x = "\nROI activity (SD)", y = "predicted rating\n") +
  plot_aes  +
  theme(legend.position = "top",
        legend.key.width=unit(2,"cm")))

H4

Do the manipulations increase neural activity in brain regions associated with self-referential processing and mentalizing?

self-referential ROI

✅ H4a: Self-focused intervention (compared to control) will increase brain activity in ROIs related to self-referential processes.

mod_h4a = lmer(self_referential ~ cond + (1 + cond | pID),
              data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h4a = table_model(mod_h4a)

table_h4a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 0.22 [0.12, 0.33] 84.20 4.33 < .001
other 0.12 [0.04, 0.19] 84.73 3.00 .004
self 0.13 [0.04, 0.22] 83.61 2.93 .004

summary

summary(mod_h4a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: self_referential ~ cond + (1 + cond | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 17259.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.4926 -0.6528 -0.0091  0.6504  3.5796 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr       
##  pID      (Intercept) 0.18580  0.4310              
##           condother   0.04599  0.2144   -0.01      
##           condself    0.08464  0.2909    0.07  0.49
##  Residual             0.97455  0.9872              
## Number of obs: 6016, groups:  pID, 85
## 
## Fixed effects:
##             Estimate Std. Error       df t value  Pr(>|t|)    
## (Intercept)  0.22400    0.05169 84.20197   4.334 0.0000404 ***
## condother    0.11680    0.03891 84.72627   3.002   0.00352 ** 
## condself     0.12998    0.04437 83.61492   2.930   0.00437 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) cndthr
## condother -0.248       
## condself  -0.166  0.488

mentalizing ROI

❌ H4b: Other-focused intervention (compared to control) will increase brain activity in ROIs related to mentalizing processes.

mod_h4b = lmer(mentalizing ~ cond + (1 + cond | pID),
              data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h4b = table_model(mod_h4b)

table_h4b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 0.27 [0.15, 0.38] 84.17 4.59 < .001
other 0.01 [-0.06, 0.09] 83.40 0.29 .774
self 0.05 [-0.03, 0.14] 84.19 1.28 .206

summary

summary(mod_h4b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: mentalizing ~ cond + (1 + cond | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 17308.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.4478 -0.6537  0.0009  0.6652  4.0095 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr       
##  pID      (Intercept) 0.24623  0.4962              
##           condother   0.03954  0.1988   -0.45      
##           condself    0.07191  0.2682   -0.25  0.87
##  Residual             0.98524  0.9926              
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##             Estimate Std. Error       df t value  Pr(>|t|)    
## (Intercept)  0.26732    0.05821 84.16977   4.593 0.0000152 ***
## condother    0.01097    0.03804 83.39828   0.288     0.774    
## condself     0.05453    0.04277 84.19090   1.275     0.206    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) cndthr
## condother -0.456       
## condself  -0.355  0.637

combined plot

predicted_h4 = ggeffects::ggpredict(mod_h4a, c("cond")) %>%
  data.frame() %>%
  mutate(atlas = "self-referential") %>%
  bind_rows(ggeffects::ggpredict(mod_h4b, c("cond")) %>%
              data.frame() %>%
              mutate(atlas = "mentalizing")) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")),
         atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

ind_data_h4 = merged %>%
  mutate(atlas = recode(atlas, "pnas_self" = "self-referential",
                        "pnas_mentalizing_nopc" = "mentalizing")) %>%
  filter(atlas %in% c("self-referential", "mentalizing")) %>%
  select(pID, cond, run, trial, atlas, parameter_estimate_std) %>%
  unique() %>%
  rename("x" = cond,
         "predicted" = parameter_estimate_std) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")),
         atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

(plot_h4 = predicted_h4 %>%
  ggplot(aes(x = x, y = predicted)) +
  stat_summary(data = ind_data_h4, aes(group = pID), fun = "mean", geom = "line",
               size = .1, color = "grey50") +
  stat_summary(aes(group = group), fun = "mean", geom = "line", size = 1) +
  geom_pointrange(aes(color = x, ymin = conf.low, ymax = conf.high), size = .75) +
  facet_grid(~atlas) +
  scale_color_manual(name = "", values = palette_condition, guide = "none") +
  scale_alpha_manual(name = "", values = c(1, .5)) +
  labs(x = "", y = "ROI activity (SD)\n") +
  plot_aes +
  theme(legend.position = c(.85, .15)))

H6

Is ROI activity positively related to sharing intentions?

self-referential ROI

✅ Stronger activity in the self-referential ROI will be related to higher sharing intentions.

mod_h6a = lmer(value ~ self_referential + (1 + self_referential | pID),
              data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h6a = table_model(mod_h6a)

table_h6a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.60 [2.52, 2.68] 85.06 67.82 < .001
self-referential 0.06 [0.03, 0.09] 83.69 4.38 < .001

summary

summary(mod_h6a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ self_referential + (1 + self_referential | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16643.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5489 -0.7231  0.1185  0.7357  2.1993 
## 
## Random effects:
##  Groups   Name             Variance Std.Dev. Corr 
##  pID      (Intercept)      0.109504 0.33091       
##           self_referential 0.004234 0.06507  -0.37
##  Residual                  0.931444 0.96511       
## Number of obs: 5937, groups:  pID, 85
## 
## Fixed effects:
##                  Estimate Std. Error       df t value             Pr(>|t|)    
## (Intercept)       2.59981    0.03833 85.06107   67.82 < 0.0000000000000002 ***
## self_referential  0.06288    0.01436 83.69285    4.38            0.0000341 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## self_rfrntl -0.260

mentalizing ROI

✅ Stronger activation in the mentalizing ROI will be related to higher sharing intentions.

mod_h6b = lmer(value ~ mentalizing + (1 + mentalizing | pID),
              data = merged_wide_alt,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h6b = table_model(mod_h6b)

table_h6b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.61 [2.53, 2.68] 84.77 68.76 < .001
mentalizing 0.04 [0.02, 0.07] 83.41 3.39 .001

summary

summary(mod_h6b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ mentalizing + (1 + mentalizing | pID)
##    Data: merged_wide_alt
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16675.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.6262 -0.7154  0.1161  0.7347  2.1016 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr
##  pID      (Intercept) 0.107434 0.32777      
##           mentalizing 0.001043 0.03229  0.08
##  Residual             0.936877 0.96792      
## Number of obs: 5941, groups:  pID, 85
## 
## Fixed effects:
##             Estimate Std. Error       df t value             Pr(>|t|)    
## (Intercept)  2.60723    0.03792 84.76823  68.760 < 0.0000000000000002 ***
## mentalizing  0.04411    0.01300 83.40864   3.394              0.00106 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## mentalizing -0.072

combined plot

vals = seq(-4.5,4.5,.1)

predicted_h6 = ggeffects::ggpredict(mod_h6a, c("self_referential [vals]")) %>%
  data.frame() %>%
  mutate(roi = "self-referential") %>%
  bind_rows(ggeffects::ggpredict(mod_h6b, c("mentalizing [vals]")) %>%
              data.frame() %>%
              mutate(roi = "mentalizing")) %>%
  mutate(roi = factor(roi, levels = c("self-referential", "mentalizing")))

ind_data_h6 = merged %>%
  select(pID, cond, run, trial, atlas, parameter_estimate_std, value) %>%
  rename("x" = parameter_estimate_std,
         "predicted" = value,
         "roi" = atlas) %>%
  mutate(roi = factor(roi, levels = c("self-referential", "mentalizing")))

predicted_h6 %>%
  ggplot(aes(x = x, y = predicted, color = roi, fill = roi)) +
  stat_smooth(data = ind_data_h6, aes(group = pID), geom ='line', method = "lm", alpha = .1, size = 1, se = FALSE) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2, color = NA) +
  geom_line(size = 2) +
  facet_grid(~roi) +
  scale_color_manual(name = "", values = palette_roi) +
  scale_fill_manual(name = "", values = palette_roi) +
  labs(y = "predicted sharing intention\n", x = "\nROI activity (SD)") +
  plot_aes +
  theme(legend.position = "none")

H7

Is there an indirect effect of the condition on sharing intentions through activity in self-referential and mentalizing ROIs?

prep data

# source functions
source("indirectMLM.R")

# create self condition dataframe
data_med_self = merged %>%
  filter(!cond == "other" & atlas == "self-referential") %>%
  mutate(cond = ifelse(cond == "self", 1, 0)) %>%
  select(pID, site, trial, cond, value, parameter_estimate) %>%
  data.frame()

# create social condition dataframe
data_med_other = merged %>%
  filter(!cond == "self" & atlas == "mentalizing") %>%
  mutate(cond = ifelse(cond == "other", 1, 0)) %>%
  select(pID, site, trial, cond, value, parameter_estimate) %>%
  data.frame()

# define variables
y_var = "value"
m_var = "parameter_estimate"

self condition

✅ H7a: The effect of self-focused intervention on sharing intention will be mediated by increased activity in the self-referential ROI.

model_name = "mediation_self"
data = data_med_self

if (file.exists(sprintf("models/model_%s_alternative.RDS", model_name))) {
  assign(get("model_name"), readRDS(sprintf("models/model_%s_alternative.RDS", model_name)))
} else {
  assign(get("model_name"), boot(data = data, statistic = indirect.mlm, R = 500,
                                 y = y_var, x = "cond", mediator = m_var, group.id = "pID",
                                 between.m = F, uncentered.x = F))
  saveRDS(eval(parse(text = model_name)), sprintf("models/model_%s_alternative.RDS", model_name))
}

indirect.mlm.summary(get(model_name))
## #### Population Covariance ####
## Covariance of Random Slopes a and b: 0.001 [-0.003, 0.01]
## 
## 
## #### Indirect Effects ####
## # Within-subject Effects
## Unbiased Estimate of Within-subjects Indirect Effect: 0.006 [0.001, 0.017]
## Biased Estimate of Within-subjects Indirect Effect: 0.005 [0, 0.012]
## Bias in Within-subjects Indirect Effect: 0.001 [0, 0.01]
## 
## 
## #### Total Effect ####
## Unbiased Estimate of Total Effect: -0.047 [-0.106, 0.008]
## Biased Total Effect of X on Y (c path): -0.044 [-0.103, 0.011]
## Bias in Total Effect: 0.002 [0, 0.007]
## 
## 
## #### Direct Effects ####
## Direct Effect of Predictor on Dependent Variable (c' path): -0.053 [-0.118, 0]
## Within-subjects Effect of Predictor on Mediator (a path for group-mean centered predictor): 0.03 [-0.001, 0.059]
## Within-subjects Effect of Mediator on Dependent Variable (b path for group-mean centered mediator): 0.173 [0.117, 0.256]

other condition

❌ H7b: The effect of other-focused intervention on sharing intention will be mediated by increased activity in the mentalizing ROI.

model_name = "mediation_other"
data = data_med_other

if (file.exists(sprintf("models/model_%s_alternative.RDS", model_name))) {
  assign(get("model_name"), readRDS(sprintf("models/model_%s_alternative.RDS", model_name)))
} else {
  assign(get("model_name"), boot(data = data, statistic = indirect.mlm, R = 500,
                                 y = y_var, x = "cond", mediator = m_var, group.id = "pID",
                                 between.m = F, uncentered.x = F))
  saveRDS(eval(parse(text = model_name)), sprintf("models/model_%s_alternative.RDS", model_name))
}

indirect.mlm.summary(get(model_name))
## #### Population Covariance ####
## Covariance of Random Slopes a and b: 0 [-0.005, 0.006]
## 
## 
## #### Indirect Effects ####
## # Within-subject Effects
## Unbiased Estimate of Within-subjects Indirect Effect: 0.003 [-0.003, 0.012]
## Biased Estimate of Within-subjects Indirect Effect: 0.003 [-0.002, 0.01]
## Bias in Within-subjects Indirect Effect: 0 [0, 0.007]
## 
## 
## #### Total Effect ####
## Unbiased Estimate of Total Effect: -0.031 [-0.093, 0.026]
## Biased Total Effect of X on Y (c path): -0.032 [-0.094, 0.029]
## Bias in Total Effect: 0.001 [0, 0.005]
## 
## 
## #### Direct Effects ####
## Direct Effect of Predictor on Dependent Variable (c' path): -0.034 [-0.099, 0.021]
## Within-subjects Effect of Predictor on Mediator (a path for group-mean centered predictor): 0.018 [-0.006, 0.039]
## Within-subjects Effect of Mediator on Dependent Variable (b path for group-mean centered mediator): 0.193 [0.148, 0.319]



combined table

table_h1a %>% mutate(DV = "H1a: Self-relevance") %>%
  bind_rows(table_h1b %>% mutate(DV = "H1b: Social relevance")) %>%
  bind_rows(table_h4a %>% mutate(DV = "H4a: Self-referential ROI")) %>%
  bind_rows(table_h4b %>% mutate(DV = "H4b: Mentalizing ROI")) %>%
  bind_rows(table_h6a %>% mutate(DV = "H6a: Sharing intention")) %>%
  bind_rows(table_h6b %>% mutate(DV = "H6b: Sharing intention")) %>%
  select(DV, everything()) %>%
  kable() %>%
  kable_styling()
DV term b [95% CI] df t p
H1a: Self-relevance intercept 2.56 [2.48, 2.63] 84.66 65.05 < .001
H1a: Self-relevance self-referential 0.03 [0.01, 0.06] 83.65 2.44 .017
H1b: Social relevance intercept 2.66 [2.58, 2.75] 84.14 64.14 < .001
H1b: Social relevance mentalizing 0.04 [0.01, 0.06] 82.92 3.10 .003
H4a: Self-referential ROI intercept 0.22 [0.12, 0.33] 84.20 4.33 < .001
H4a: Self-referential ROI other 0.12 [0.04, 0.19] 84.73 3.00 .004
H4a: Self-referential ROI self 0.13 [0.04, 0.22] 83.61 2.93 .004
H4b: Mentalizing ROI intercept 0.27 [0.15, 0.38] 84.17 4.59 < .001
H4b: Mentalizing ROI other 0.01 [-0.06, 0.09] 83.40 0.29 .774
H4b: Mentalizing ROI self 0.05 [-0.03, 0.14] 84.19 1.28 .206
H6a: Sharing intention intercept 2.60 [2.52, 2.68] 85.06 67.82 < .001
H6a: Sharing intention self-referential 0.06 [0.03, 0.09] 83.69 4.38 < .001
H6b: Sharing intention intercept 2.61 [2.53, 2.68] 84.77 68.76 < .001
H6b: Sharing intention mentalizing 0.04 [0.02, 0.07] 83.41 3.39 .001

moderation by article topic

These analyses explore whether the analyses reported in study 1 of the main manuscript are moderated by article topic (health or climate).

H1

Are the relationships between ROI activity and self and social relevance ratings moderated by article topic?

self-referential ROI

There is a main effect of topic, such that health articles elicited greater activity in the self-referential ROI compared to climate articles.

These data are not consistent with moderation by topic.

mod_h1a =  lmer(self_relevance ~ self_referential * topic + (1 | pID),
               data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h1a = table_model(mod_h1a)

table_h1a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.49 [2.41, 2.57] 102.39 61.76 < .001
self-referential 0.04 [0.00, 0.07] 6015.96 2.11 .035
topic (health) 0.14 [0.09, 0.19] 5933.30 5.56 < .001
self-referential x topic (health) 0.02 [-0.02, 0.06] 5948.49 0.86 .390

summary

summary(mod_h1a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: self_relevance ~ self_referential * topic + (1 | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16752.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5459 -0.6971  0.1408  0.6762  2.3717 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  pID      (Intercept) 0.1124   0.3352  
##  Residual             0.9136   0.9558  
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                                Estimate Std. Error         df t value
## (Intercept)                     2.49156    0.04034  102.39475  61.758
## self_referential                0.03518    0.01665 6015.95934   2.113
## topichealth                     0.13827    0.02486 5933.29867   5.561
## self_referential:topichealth    0.01925    0.02238 5948.48840   0.860
##                                          Pr(>|t|)    
## (Intercept)                  < 0.0000000000000002 ***
## self_referential                           0.0346 *  
## topichealth                          0.0000000279 ***
## self_referential:topichealth               0.3895    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) slf_rf tpchlt
## self_rfrntl -0.034              
## topichealth -0.304  0.043       
## slf_rfrntl:  0.023 -0.676 -0.122

mentalizing ROI

There is a main effect of topic, such that health articles elicited greater activity in the mentalizing ROI compared to climate articles.

These data are not consistent with moderation by topic.

mod_h1b = lmer(social_relevance ~ mentalizing * topic + (1 + mentalizing | pID),
               data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h1b = table_model(mod_h1b)

table_h1b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.51 [2.43, 2.60] 98.00 58.11 < .001
mentalizing 0.04 [0.01, 0.07] 238.75 2.37 .019
topic (health) 0.29 [0.25, 0.34] 5928.14 12.30 < .001
mentalizing x topic (health) 0.01 [-0.03, 0.05] 5892.58 0.55 .586

summary

summary(mod_h1b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: social_relevance ~ mentalizing * topic + (1 + mentalizing | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 15692.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7751 -0.7054  0.1247  0.6758  2.7757 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr 
##  pID      (Intercept) 0.134968 0.36738       
##           mentalizing 0.001087 0.03297  -0.20
##  Residual             0.761555 0.87267       
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                           Estimate Std. Error         df t value
## (Intercept)                2.51144    0.04322   97.99658  58.112
## mentalizing                0.03671    0.01550  238.74662   2.369
## topichealth                0.29267    0.02380 5928.13941  12.299
## mentalizing:topichealth    0.01119    0.02054 5892.57757   0.545
##                                    Pr(>|t|)    
## (Intercept)             <0.0000000000000002 ***
## mentalizing                          0.0186 *  
## topichealth             <0.0000000000000002 ***
## mentalizing:topichealth              0.5857    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) mntlzn tpchlt
## mentalizing -0.152              
## topichealth -0.267  0.175       
## mntlzng:tpc  0.075 -0.649 -0.321

combined plot

predicted = ggeffects::ggpredict(mod_h1a, c("self_referential [-4.5:5]", "topic")) %>%
  data.frame() %>%
  mutate(roi = "self-referential",
         variable = "self-relevance") %>%
  bind_rows(ggeffects::ggpredict(mod_h1b, c("mentalizing [-4.5:5]", "topic")) %>%
              data.frame() %>%
              mutate(roi = "mentalizing",
                     variable = "social relevance"))

ind_data = merged_wide %>%
  select(topic, pID, trial, contains("relevance"), mentalizing, self_referential) %>%
  rename("self-referential" = self_referential,
         "group" = topic) %>%
  gather(variable, predicted, contains("relevance")) %>%
  mutate(variable = gsub("self_relevance", "self-relevance", variable),
         variable = gsub("social_relevance", "social relevance", variable)) %>%
  gather(roi, x, mentalizing, `self-referential`) %>%
  filter(!(variable == "self-relevance" & roi == "mentalizing") & ! (variable == "social relevance" & roi == "self-referential"))

(plot_h1 = predicted %>%
  ggplot(aes(x, predicted, color = group, fill = group)) +
  stat_smooth(data = ind_data, aes(group = interaction(pID, group)), geom ='line', method = "lm", alpha = .1, size = 1, se = FALSE) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .3, color = NA) +
  geom_line(size = 2) +
  facet_grid(~variable) +
  scale_color_manual(name = "", values = palette_topic) +
  scale_fill_manual(name = "", values = palette_topic) +
  labs(x = "\nROI activity (SD)", y = "predicted rating\n") +
  plot_aes  +
  theme(legend.position = "top",
        legend.key.width=unit(2,"cm")))

H2

Are the effects of the experimental manipulations on relevance moderated by article topic?

self-relevance

There is a main effect of topic such that health articles are rated as more self-relevant than climate articles.

The was also an interaction such that the effect of the self-focused condition on self-relevance was weaker for health articles.

mod_h2a = lmer(self_relevance ~ cond * topic + (1 | pID),
               data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h2a = table_model(mod_h2a)

table_h2a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.44 [2.35, 2.53] 192.41 51.76 < .001
other 0.04 [-0.04, 0.13] 5930.30 1.04 .299
self 0.12 [0.04, 0.20] 5930.32 2.80 .005
topic (health) 0.23 [0.14, 0.31] 5930.31 5.35 < .001
other x topic (health) -0.08 [-0.19, 0.04] 5930.38 -1.26 .208
self x topic (health) -0.17 [-0.29, -0.05] 5930.36 -2.86 .004

summary

summary(mod_h2a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: self_relevance ~ cond * topic + (1 | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16762.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5028 -0.6952  0.1359  0.6721  2.4245 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  pID      (Intercept) 0.1116   0.3341  
##  Residual             0.9147   0.9564  
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                         Estimate Std. Error         df t value
## (Intercept)              2.44002    0.04714  192.41193  51.764
## condother                0.04433    0.04269 5930.29698   1.038
## condself                 0.11933    0.04268 5930.32430   2.796
## topichealth              0.22857    0.04269 5930.30913   5.354
## condother:topichealth   -0.07600    0.06038 5930.38150  -1.259
## condself:topichealth    -0.17266    0.06040 5930.35977  -2.859
##                                   Pr(>|t|)    
## (Intercept)           < 0.0000000000000002 ***
## condother                          0.29915    
## condself                           0.00519 ** 
## topichealth                   0.0000000892 ***
## condother:topichealth              0.20821    
## condself:topichealth               0.00427 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) cndthr cndslf tpchlt cndth:
## condother   -0.451                            
## condself    -0.452  0.499                     
## topichealth -0.451  0.499  0.499              
## cndthr:tpch  0.319 -0.707 -0.353 -0.707       
## cndslf:tpch  0.319 -0.352 -0.707 -0.707  0.500

social relevance

There is a main effect of topic such that health articles are rated as more socially relevant than climate articles.

These data are not consistent with moderation by topic.

mod_h2b = lmer(social_relevance ~ cond * topic + (1 | pID),
               data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h2b = table_model(mod_h2b)

table_h2b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.49 [2.39, 2.58] 158.32 51.66 < .001
other 0.03 [-0.05, 0.11] 5930.24 0.79 .431
self 0.07 [-0.00, 0.15] 5930.26 1.91 .056
topic (health) 0.31 [0.23, 0.39] 5930.24 7.97 < .001
other x topic (health) 0.03 [-0.08, 0.14] 5930.30 0.52 .606
self x topic (health) -0.06 [-0.17, 0.05] 5930.28 -1.06 .289

summary

summary(mod_h2b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: social_relevance ~ cond * topic + (1 | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 15706.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7954 -0.7085  0.1160  0.6672  2.7648 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  pID      (Intercept) 0.1328   0.3644  
##  Residual             0.7640   0.8741  
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                         Estimate Std. Error         df t value
## (Intercept)              2.48871    0.04818  158.32368  51.655
## condother                0.03072    0.03902 5930.23528   0.787
## condself                 0.07450    0.03901 5930.25561   1.910
## topichealth              0.31089    0.03902 5930.24394   7.968
## condother:topichealth    0.02844    0.05518 5930.29664   0.515
## condself:topichealth    -0.05854    0.05520 5930.28121  -1.061
##                                   Pr(>|t|)    
## (Intercept)           < 0.0000000000000002 ***
## condother                           0.4311    
## condself                            0.0562 .  
## topichealth            0.00000000000000191 ***
## condother:topichealth               0.6063    
## condself:topichealth                0.2889    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) cndthr cndslf tpchlt cndth:
## condother   -0.404                            
## condself    -0.404  0.499                     
## topichealth -0.404  0.499  0.499              
## cndthr:tpch  0.285 -0.707 -0.353 -0.707       
## cndslf:tpch  0.285 -0.352 -0.707 -0.707  0.500

combined plot

predicted_h2 = ggeffects::ggpredict(mod_h2a, c("cond", "topic")) %>%
  data.frame() %>%
  mutate(model = "self-relevance") %>%
  bind_rows(ggeffects::ggpredict(mod_h2b, c("cond", "topic")) %>%
              data.frame() %>%
              mutate(model = "social relevance")) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")))

ind_data_h2 = merged_wide %>%
  rename("x" = cond,
         "group" = topic) %>%
  gather(model, predicted, self_relevance, social_relevance) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")),
         model = gsub("self_relevance", "self-relevance", model),
         model = gsub("social_relevance", "social relevance", model))
  
(plot_h2 = predicted_h2 %>%
  ggplot(aes(x = x, y = predicted, color = group)) +
  stat_summary(data = ind_data_h2, aes(group = interaction(pID, group)), fun = "mean", geom = "line", size = .1) +
  stat_summary(aes(group = group), fun = "mean", geom = "line", size = 1, position = position_dodge(.1)) +
  geom_pointrange(aes(ymin = conf.low, ymax = conf.high, group = group),
                  size = .75, position = position_dodge(.1)) +
  facet_grid(~model) +
  scale_color_manual(name = "", values = palette_topic) +
  labs(x = "", y = "predicted rating\n") +
  plot_aes +
  theme(legend.position = c(.85, .15)))

H3

Are the relationships between self and social relevance and sharing intentions moderated by article topic?

The relationship between self-relevance and sharing intentions was not moderated by topic.

However, the relationship between social relevance and sharing intentions was slightly stronger for health articles compared to climate articles.

mod_h3 = lmer(value ~ self_relevance * topic + social_relevance * topic + (1 + self_relevance + social_relevance | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

plot

predicted = ggeffects::ggpredict(mod_h3, c("self_relevance", "topic")) %>%
  data.frame() %>%
  mutate(variable = "self-relevance") %>%
  bind_rows(ggeffects::ggpredict(mod_h3, c("social_relevance", "topic")) %>%
              data.frame() %>%
              mutate(variable = "social relevance"))

points = merged_wide %>%
  rename("self-referential" = self_referential,
         "predicted" = value,
         "group" = topic) %>%
  gather(variable, x, contains("relevance")) %>%
  mutate(variable = gsub("self_relevance", "self-relevance", variable),
         variable = gsub("social_relevance", "social relevance", variable))

(plot_rel_sharing = predicted %>%
  ggplot(aes(x, predicted, color = group, fill = group)) +
  stat_smooth(data = points, aes(group = interaction(pID, group)), geom ='line', method = "lm", alpha = .1, size = 1, se = FALSE) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2, color = NA) +
  geom_line(size = 2) +
  facet_grid(~variable) +
  scale_color_manual(name = "", values = palette_topic) +
  scale_fill_manual(name = "", values = palette_topic) +
  labs(x = "\nrating", y = "predicted sharing intention\n") +
  plot_aes)

model table

table_h3 = table_model(mod_h3)

table_h3 %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 1.26 [1.11, 1.40] 116.06 17.66 < .001
self-relevance 0.29 [0.24, 0.34] 275.53 10.79 < .001
topic (health) -0.11 [-0.24, 0.03] 5860.61 -1.51 .130
social relevance 0.21 [0.14, 0.27] 179.99 6.54 < .001
self-relevance x topic (health) 0.03 [-0.03, 0.09] 5557.57 1.04 .300
topic (health) x social relevance 0.06 [0.00, 0.13] 5632.89 1.98 .048

summary

summary(mod_h3)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ self_relevance * topic + social_relevance * topic + (1 +  
##     self_relevance + social_relevance | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 14876.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.2184 -0.6968  0.0616  0.6860  3.0642 
## 
## Random effects:
##  Groups   Name             Variance Std.Dev. Corr       
##  pID      (Intercept)      0.22603  0.4754              
##           self_relevance   0.01285  0.1134   -0.27      
##           social_relevance 0.03180  0.1783   -0.58 -0.53
##  Residual                  0.67979  0.8245              
## Number of obs: 5941, groups:  pID, 85
## 
## Fixed effects:
##                                Estimate Std. Error         df t value
## (Intercept)                     1.25511    0.07109  116.06407  17.656
## self_relevance                  0.29150    0.02702  275.53284  10.787
## topichealth                    -0.10627    0.07022 5860.60678  -1.513
## social_relevance                0.20677    0.03163  179.99317   6.537
## self_relevance:topichealth      0.03147    0.03034 5557.57118   1.037
## topichealth:social_relevance    0.06340    0.03201 5632.89424   1.981
##                                          Pr(>|t|)    
## (Intercept)                  < 0.0000000000000002 ***
## self_relevance               < 0.0000000000000002 ***
## topichealth                                0.1303    
## social_relevance                   0.000000000627 ***
## self_relevance:topichealth                 0.2997    
## topichealth:social_relevance               0.0476 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) slf_rl tpchlt scl_rl slf_r:
## self_relvnc -0.268                            
## topichealth -0.438  0.182                     
## socil_rlvnc -0.482 -0.657  0.192              
## slf_rlvnc:t  0.154 -0.679 -0.293  0.435       
## tpchlth:sc_  0.196  0.477 -0.467 -0.575 -0.675

H4

Are the effects of the experimental manipulations on ROI activity moderated by article topic?

self-referential ROI

There is a main effect of topic, such that health articles elicited greater activity in the self-referential ROI compared to climate articles.

These data are not consistent with moderation by topic.

mod_h4a = lmer(self_referential ~ cond * topic + (1 + cond | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h4a = table_model(mod_h4a)

table_h4a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 0.01 [-0.11, 0.13] 111.47 0.12 .909
other 0.09 [-0.00, 0.19] 229.04 1.89 .060
self 0.13 [0.03, 0.24] 194.36 2.51 .013
topic (health) 0.14 [0.05, 0.23] 5765.56 3.15 .002
other x topic (health) -0.01 [-0.13, 0.11] 5765.13 -0.17 .867
self x topic (health) -0.08 [-0.20, 0.04] 5765.23 -1.30 .192

summary

summary(mod_h4a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: self_referential ~ cond * topic + (1 + cond | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 17294.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.8462 -0.6484  0.0071  0.6462  3.5919 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr       
##  pID      (Intercept) 0.22957  0.4791              
##           condother   0.04296  0.2073   -0.16      
##           condself    0.07414  0.2723   -0.06  0.58
##  Residual             0.97687  0.9884              
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                          Estimate  Std. Error          df t value Pr(>|t|)   
## (Intercept)              0.006975    0.060591  111.469116   0.115  0.90856   
## condother                0.093636    0.049521  229.036563   1.891  0.05991 . 
## condself                 0.133495    0.053089  194.356460   2.515  0.01273 * 
## topichealth              0.138852    0.044121 5765.562012   3.147  0.00166 **
## condother:topichealth   -0.010460    0.062405 5765.129221  -0.168  0.86690   
## condself:topichealth    -0.081383    0.062421 5765.229257  -1.304  0.19236   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) cndthr cndslf tpchlt cndth:
## condother   -0.384                            
## condself    -0.331  0.516                     
## topichealth -0.363  0.444  0.414              
## cndthr:tpch  0.257 -0.630 -0.293 -0.707       
## cndslf:tpch  0.257 -0.314 -0.587 -0.707  0.500

mentalizing ROI

There is a main effect of topic, such that health articles elicited greater activity in the mentalizing ROI compared to climate articles.

These data are not consistent with moderation by topic.

mod_h4b = lmer(mentalizing ~ cond * topic + (1 + cond | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h4b = table_model(mod_h4b)

table_h4b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 0.26 [0.15, 0.38] 113.11 4.46 < .001
other 0.06 [-0.03, 0.16] 240.17 1.28 .202
self 0.11 [0.01, 0.21] 199.26 2.08 .038
topic (health) 0.11 [0.03, 0.20] 5765.52 2.59 .010
other x topic (health) 0.00 [-0.12, 0.12] 5765.10 0.04 .972
self x topic (health) -0.06 [-0.18, 0.06] 5765.16 -0.96 .339

summary

summary(mod_h4b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: mentalizing ~ cond * topic + (1 + cond | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 17303.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.7051 -0.6544  0.0174  0.6623  3.2816 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr       
##  pID      (Intercept) 0.21622  0.4650              
##           condother   0.03545  0.1883   -0.15      
##           condself    0.06887  0.2624   -0.03  0.61
##  Residual             0.98059  0.9902              
## Number of obs: 6020, groups:  pID, 85
## 
## Fixed effects:
##                          Estimate  Std. Error          df t value  Pr(>|t|)    
## (Intercept)              0.264690    0.059311  113.112427   4.463 0.0000192 ***
## condother                0.062299    0.048695  240.165739   1.279    0.2020    
## condself                 0.109549    0.052572  199.264152   2.084    0.0385 *  
## topichealth              0.114532    0.044205 5765.519256   2.591    0.0096 ** 
## condother:topichealth    0.002201    0.062523 5765.102401   0.035    0.9719    
## condself:topichealth    -0.059863    0.062540 5765.158731  -0.957    0.3385    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) cndthr cndslf tpchlt cndth:
## condother   -0.392                            
## condself    -0.328  0.519                     
## topichealth -0.371  0.452  0.419              
## cndthr:tpch  0.263 -0.642 -0.296 -0.707       
## cndslf:tpch  0.263 -0.320 -0.594 -0.707  0.500

combined plot

predicted_h4 = ggeffects::ggpredict(mod_h4a, c("cond", "topic")) %>%
  data.frame() %>%
  mutate(atlas = "self-referential") %>%
  bind_rows(ggeffects::ggpredict(mod_h4b, c("cond", "topic")) %>%
              data.frame() %>%
              mutate(atlas = "mentalizing")) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")),
         atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

ind_data_h4 = merged %>%
  filter(atlas %in% c("self-referential", "mentalizing")) %>%
  select(topic, pID, cond, run, trial, atlas, parameter_estimate_std) %>%
  unique() %>%
  rename("x" = cond,
         "predicted" = parameter_estimate_std,
         "group" = topic) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")),
         atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

(plot_h4 = predicted_h4 %>%
  ggplot(aes(x = x, y = predicted, color = group)) +
  stat_summary(data = ind_data_h4, aes(group = interaction(pID, group)), fun = "mean", geom = "line", size = .1) +
  stat_summary(aes(group = group), fun = "mean", geom = "line", size = 1, position = position_dodge(.1)) +
  geom_pointrange(aes(ymin = conf.low, ymax = conf.high, group = group),
                  size = .75, position = position_dodge(.1)) +
  facet_grid(~atlas) +
  scale_color_manual(name = "", values = palette_topic) +
  labs(x = "", y = "ROI activity (SD)\n") +
  plot_aes +
  theme(legend.position = c(.85, .15)))

H5

Are the effect of the experimental manipulations on sharing intentions moderated by article topic?

There is a main effect of topic, such that health articles have higher sharing intentions than climate articles.

These data are not consistent with moderation by topic.

mod_h5 = lmer(value ~ cond * topic + (1 | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

plot

predicted_h5 = ggeffects::ggpredict(mod_h5, c("cond", "topic")) %>%
  data.frame() %>%
  mutate(x = factor(x, levels = c("self", "control", "other")))

ind_data_h5 = merged_wide %>%
  rename("x" = cond,
         "predicted" = value,
         "group" = topic) %>%
  mutate(x = factor(x, levels = c("self", "control", "other")))
  
predicted_h5 %>%
  ggplot(aes(x = x, y = predicted, color = group)) +
  stat_summary(data = ind_data_h5, aes(group = interaction(pID, group)), fun = "mean", geom = "line", size = .1) +
  stat_summary(aes(group = group), fun = "mean", geom = "line", size = 1, position = position_dodge(.1)) +
  geom_pointrange(aes(ymin = conf.low, ymax = conf.high, group = group),
                  size = .75, position = position_dodge(.1)) +
  scale_color_manual(name = "", values = palette_topic) +
  labs(x = "", y = "predicted sharing intention\n") +
  plot_aes +
  theme(legend.position = c(.85, .15))

model table

table_h5 = table_model(mod_h5)

table_h5 %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.53 [2.43, 2.62] 201.73 54.11 < .001
other -0.06 [-0.14, 0.03] 5851.60 -1.29 .198
self -0.04 [-0.13, 0.04] 5851.56 -0.94 .347
topic (health) 0.24 [0.15, 0.32] 5851.55 5.52 < .001
other x topic (health) 0.05 [-0.07, 0.17] 5851.68 0.76 .450
self x topic (health) -0.01 [-0.13, 0.11] 5851.65 -0.16 .873

summary

summary(mod_h5)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ cond * topic + (1 | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16602.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5098 -0.7361  0.1117  0.7452  2.1738 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  pID      (Intercept) 0.1062   0.3259  
##  Residual             0.9244   0.9614  
## Number of obs: 5941, groups:  pID, 85
## 
## Fixed effects:
##                         Estimate Std. Error         df t value
## (Intercept)              2.52661    0.04670  201.73020  54.107
## condother               -0.05570    0.04322 5851.59600  -1.289
## condself                -0.04068    0.04324 5851.56406  -0.941
## topichealth              0.23829    0.04319 5851.54792   5.518
## condother:topichealth    0.04617    0.06109 5851.68484   0.756
## condself:topichealth    -0.00979    0.06112 5851.64700  -0.160
##                                   Pr(>|t|)    
## (Intercept)           < 0.0000000000000002 ***
## condother                            0.198    
## condself                             0.347    
## topichealth                   0.0000000358 ***
## condother:topichealth                0.450    
## condself:topichealth                 0.873    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) cndthr cndslf tpchlt cndth:
## condother   -0.461                            
## condself    -0.461  0.498                     
## topichealth -0.462  0.499  0.499              
## cndthr:tpch  0.326 -0.707 -0.352 -0.707       
## cndslf:tpch  0.326 -0.353 -0.707 -0.707  0.499

H6

Are the relationships between ROI activity positively and sharing intentions moderated by article topic?

self-referential ROI

There is a main effect of topic, such that health articles have higher sharing intentions than climate articles.

These data are not consistent with moderation by topic.

mod_h6a = lmer(value ~ self_referential * topic + (1 + self_referential | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h6a = table_model(mod_h6a)

table_h6a %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.49 [2.41, 2.57] 103.51 62.47 < .001
self-referential 0.07 [0.04, 0.11] 242.71 4.15 < .001
topic (health) 0.24 [0.19, 0.29] 5850.96 9.61 < .001
self-referential x topic (health) 0.01 [-0.04, 0.05] 5842.67 0.34 .731

summary

summary(mod_h6a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ self_referential * topic + (1 + self_referential | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16560.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5073 -0.7312  0.1088  0.7522  2.1642 
## 
## Random effects:
##  Groups   Name             Variance Std.Dev. Corr 
##  pID      (Intercept)      0.107849 0.32840       
##           self_referential 0.001754 0.04188  -0.32
##  Residual                  0.916753 0.95747       
## Number of obs: 5941, groups:  pID, 85
## 
## Fixed effects:
##                                 Estimate  Std. Error          df t value
## (Intercept)                     2.487268    0.039815  103.508777  62.470
## self_referential                0.072425    0.017469  242.708050   4.146
## topichealth                     0.241218    0.025106 5850.963104   9.608
## self_referential:topichealth    0.007782    0.022628 5842.669873   0.344
##                                          Pr(>|t|)    
## (Intercept)                  < 0.0000000000000002 ***
## self_referential                        0.0000468 ***
## topichealth                  < 0.0000000000000002 ***
## self_referential:topichealth                0.731    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) slf_rf tpchlt
## self_rfrntl -0.113              
## topichealth -0.312  0.045       
## slf_rfrntl:  0.026 -0.655 -0.126

mentalizing ROI

There is a main effect of topic, such that health articles have higher sharing intentions than climate articles.

These data are not consistent with moderation by topic.

mod_h6b = lmer(value ~ mentalizing * topic + (1 + mentalizing | pID),
              data = merged_wide,
              control = lmerControl(optimizer = "bobyqa"))

model table

table_h6b = table_model(mod_h6b)

table_h6b %>%
    kable()  %>%
    kableExtra::kable_styling()
term b [95% CI] df t p
intercept 2.47 [2.39, 2.55] 105.78 61.42 < .001
mentalizing 0.07 [0.04, 0.11] 242.79 4.29 < .001
topic (health) 0.25 [0.20, 0.30] 5846.05 9.42 < .001
mentalizing x topic (health) -0.01 [-0.06, 0.03] 5827.12 -0.50 .616

summary

summary(mod_h6b)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ mentalizing * topic + (1 + mentalizing | pID)
##    Data: merged_wide
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 16569.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.4736 -0.7340  0.1136  0.7598  2.1242 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr 
##  pID      (Intercept) 0.10830  0.32909       
##           mentalizing 0.00101  0.03178  -0.20
##  Residual             0.91878  0.95853       
## Number of obs: 5941, groups:  pID, 85
## 
## Fixed effects:
##                           Estimate Std. Error         df t value
## (Intercept)                2.47032    0.04022  105.78117  61.418
## mentalizing                0.07312    0.01705  242.78508   4.288
## topichealth                0.24818    0.02634 5846.05395   9.421
## mentalizing:topichealth   -0.01139    0.02272 5827.12098  -0.501
##                                     Pr(>|t|)    
## (Intercept)             < 0.0000000000000002 ***
## mentalizing                         0.000026 ***
## topichealth             < 0.0000000000000002 ***
## mentalizing:topichealth                0.616    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) mntlzn tpchlt
## mentalizing -0.170              
## topichealth -0.318  0.182       
## mntlzng:tpc  0.092 -0.656 -0.325

combined plot

vals = seq(-4.5,4.5,.1)

predicted_h6 = ggeffects::ggpredict(mod_h6a, c("self_referential [vals]", "topic")) %>%
  data.frame() %>%
  mutate(atlas = "self-referential") %>%
  bind_rows(ggeffects::ggpredict(mod_h6b, c("mentalizing [vals]", "topic")) %>%
              data.frame() %>%
              mutate(atlas = "mentalizing")) %>%
  mutate(atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

ind_data_h6 = merged %>%
  filter(atlas %in% c("self-referential", "mentalizing")) %>%
  select(topic, pID, cond, run, trial, atlas, parameter_estimate_std, value) %>%
  rename("x" = parameter_estimate_std,
         "predicted" = value,
         "group" = topic) %>%
  mutate(atlas = factor(atlas, levels = c("self-referential", "mentalizing")))

predicted_h6 %>%
  ggplot(aes(x = x, y = predicted, color = group, fill = group)) +
  stat_smooth(data = ind_data_h6, aes(group = interaction(pID, group)), geom ='line', method = "lm", alpha = .1, size = 1, se = FALSE) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2, color = NA) +
  geom_line(size = 2) +
  facet_grid(~atlas) +
  scale_color_manual(name = "", values = palette_topic) +
  scale_fill_manual(name = "", values = palette_topic) +
  labs(y = "predicted sharing intention\n", x = "\nROI activity (SD)") +
  plot_aes +
  theme(legend.position = "top")

combined table

table_h1a %>% mutate(DV = "H1a: Self-relevance") %>%
  bind_rows(table_h1b %>% mutate(DV = "H1b: Social relevance")) %>%
  bind_rows(table_h2a %>% mutate(DV = "H2a: Self-relevance")) %>%
  bind_rows(table_h2b %>% mutate(DV = "H2b: Social relevance")) %>%
  bind_rows(table_h3 %>% mutate(DV = "H3a-b: Sharing intention")) %>%
  bind_rows(table_h4a %>% mutate(DV = "H4a: Self-referential ROI")) %>%
  bind_rows(table_h4b %>% mutate(DV = "H4b: Mentalizing ROI")) %>%
  bind_rows(table_h5 %>% mutate(DV = "H5: Sharing intention")) %>%
  bind_rows(table_h6a %>% mutate(DV = "H6a: Sharing intention")) %>%
  bind_rows(table_h6b %>% mutate(DV = "H6b: Sharing intention")) %>%
  select(DV, everything()) %>%
  kable() %>%
  kable_styling()
DV term b [95% CI] df t p
H1a: Self-relevance intercept 2.49 [2.41, 2.57] 102.39 61.76 < .001
H1a: Self-relevance self-referential 0.04 [0.00, 0.07] 6015.96 2.11 .035
H1a: Self-relevance topic (health) 0.14 [0.09, 0.19] 5933.30 5.56 < .001
H1a: Self-relevance self-referential x topic (health) 0.02 [-0.02, 0.06] 5948.49 0.86 .390
H1b: Social relevance intercept 2.51 [2.43, 2.60] 98.00 58.11 < .001
H1b: Social relevance mentalizing 0.04 [0.01, 0.07] 238.75 2.37 .019
H1b: Social relevance topic (health) 0.29 [0.25, 0.34] 5928.14 12.30 < .001
H1b: Social relevance mentalizing x topic (health) 0.01 [-0.03, 0.05] 5892.58 0.55 .586
H2a: Self-relevance intercept 2.44 [2.35, 2.53] 192.41 51.76 < .001
H2a: Self-relevance other 0.04 [-0.04, 0.13] 5930.30 1.04 .299
H2a: Self-relevance self 0.12 [0.04, 0.20] 5930.32 2.80 .005
H2a: Self-relevance topic (health) 0.23 [0.14, 0.31] 5930.31 5.35 < .001
H2a: Self-relevance other x topic (health) -0.08 [-0.19, 0.04] 5930.38 -1.26 .208
H2a: Self-relevance self x topic (health) -0.17 [-0.29, -0.05] 5930.36 -2.86 .004
H2b: Social relevance intercept 2.49 [2.39, 2.58] 158.32 51.66 < .001
H2b: Social relevance other 0.03 [-0.05, 0.11] 5930.24 0.79 .431
H2b: Social relevance self 0.07 [-0.00, 0.15] 5930.26 1.91 .056
H2b: Social relevance topic (health) 0.31 [0.23, 0.39] 5930.24 7.97 < .001
H2b: Social relevance other x topic (health) 0.03 [-0.08, 0.14] 5930.30 0.52 .606
H2b: Social relevance self x topic (health) -0.06 [-0.17, 0.05] 5930.28 -1.06 .289
H3a-b: Sharing intention intercept 1.26 [1.11, 1.40] 116.06 17.66 < .001
H3a-b: Sharing intention self-relevance 0.29 [0.24, 0.34] 275.53 10.79 < .001
H3a-b: Sharing intention topic (health) -0.11 [-0.24, 0.03] 5860.61 -1.51 .130
H3a-b: Sharing intention social relevance 0.21 [0.14, 0.27] 179.99 6.54 < .001
H3a-b: Sharing intention self-relevance x topic (health) 0.03 [-0.03, 0.09] 5557.57 1.04 .300
H3a-b: Sharing intention topic (health) x social relevance 0.06 [0.00, 0.13] 5632.89 1.98 .048
H4a: Self-referential ROI intercept 0.01 [-0.11, 0.13] 111.47 0.12 .909
H4a: Self-referential ROI other 0.09 [-0.00, 0.19] 229.04 1.89 .060
H4a: Self-referential ROI self 0.13 [0.03, 0.24] 194.36 2.51 .013
H4a: Self-referential ROI topic (health) 0.14 [0.05, 0.23] 5765.56 3.15 .002
H4a: Self-referential ROI other x topic (health) -0.01 [-0.13, 0.11] 5765.13 -0.17 .867
H4a: Self-referential ROI self x topic (health) -0.08 [-0.20, 0.04] 5765.23 -1.30 .192
H4b: Mentalizing ROI intercept 0.26 [0.15, 0.38] 113.11 4.46 < .001
H4b: Mentalizing ROI other 0.06 [-0.03, 0.16] 240.17 1.28 .202
H4b: Mentalizing ROI self 0.11 [0.01, 0.21] 199.26 2.08 .038
H4b: Mentalizing ROI topic (health) 0.11 [0.03, 0.20] 5765.52 2.59 .010
H4b: Mentalizing ROI other x topic (health) 0.00 [-0.12, 0.12] 5765.10 0.04 .972
H4b: Mentalizing ROI self x topic (health) -0.06 [-0.18, 0.06] 5765.16 -0.96 .339
H5: Sharing intention intercept 2.53 [2.43, 2.62] 201.73 54.11 < .001
H5: Sharing intention other -0.06 [-0.14, 0.03] 5851.60 -1.29 .198
H5: Sharing intention self -0.04 [-0.13, 0.04] 5851.56 -0.94 .347
H5: Sharing intention topic (health) 0.24 [0.15, 0.32] 5851.55 5.52 < .001
H5: Sharing intention other x topic (health) 0.05 [-0.07, 0.17] 5851.68 0.76 .450
H5: Sharing intention self x topic (health) -0.01 [-0.13, 0.11] 5851.65 -0.16 .873
H6a: Sharing intention intercept 2.49 [2.41, 2.57] 103.51 62.47 < .001
H6a: Sharing intention self-referential 0.07 [0.04, 0.11] 242.71 4.15 < .001
H6a: Sharing intention topic (health) 0.24 [0.19, 0.29] 5850.96 9.61 < .001
H6a: Sharing intention self-referential x topic (health) 0.01 [-0.04, 0.05] 5842.67 0.34 .731
H6b: Sharing intention intercept 2.47 [2.39, 2.55] 105.78 61.42 < .001
H6b: Sharing intention mentalizing 0.07 [0.04, 0.11] 242.79 4.29 < .001
H6b: Sharing intention topic (health) 0.25 [0.20, 0.30] 5846.05 9.42 < .001
H6b: Sharing intention mentalizing x topic (health) -0.01 [-0.06, 0.03] 5827.12 -0.50 .616

cite packages

report::cite_packages()
##   - Angelo Canty and Brian Ripley (2021). boot: Bootstrap R (S-Plus) Functions. R package version 1.3-28.
##   - Douglas Bates and Martin Maechler (2021). Matrix: Sparse and Dense Matrix Classes and Methods. R package version 1.3-4. https://CRAN.R-project.org/package=Matrix
##   - Douglas Bates, Martin Maechler, Ben Bolker, Steve Walker (2015). Fitting Linear Mixed-Effects Models Using lme4. Journal of Statistical Software, 67(1), 1-48. doi:10.18637/jss.v067.i01.
##   - H. Wickham. ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York, 2016.
##   - Hadley Wickham (2019). stringr: Simple, Consistent Wrappers for Common String Operations. R package version 1.4.0. https://CRAN.R-project.org/package=stringr
##   - Hadley Wickham (2021). forcats: Tools for Working with Categorical Variables (Factors). R package version 0.5.1. https://CRAN.R-project.org/package=forcats
##   - Hadley Wickham and Maximilian Girlich (2022). tidyr: Tidy Messy Data. R package version 1.2.0. https://CRAN.R-project.org/package=tidyr
##   - Hadley Wickham, Jennifer Bryan and Malcolm Barrett (2021). usethis: Automate Package and Project Setup. R package version 2.1.5. https://CRAN.R-project.org/package=usethis
##   - Hadley Wickham, Jim Hester and Jennifer Bryan (2022). readr: Read Rectangular Text Data. R package version 2.1.2. https://CRAN.R-project.org/package=readr
##   - Hadley Wickham, Jim Hester, Winston Chang and Jennifer Bryan (2021). devtools: Tools to Make Developing R Packages Easier. R package version 2.4.3. https://CRAN.R-project.org/package=devtools
##   - Hadley Wickham, Romain François, Lionel Henry and Kirill Müller (2022). dplyr: A Grammar of Data Manipulation. R package version 1.0.9. https://CRAN.R-project.org/package=dplyr
##   - Hao Zhu (2021). kableExtra: Construct Complex Table with 'kable' and Pipe Syntax. R package version 1.3.4. https://CRAN.R-project.org/package=kableExtra
##   - Jim Hester, Hadley Wickham and Gábor Csárdi (2021). fs: Cross-Platform File System Operations Based on 'libuv'. R package version 1.5.2. https://CRAN.R-project.org/package=fs
##   - Kirill Müller and Hadley Wickham (2022). tibble: Simple Data Frames. R package version 3.1.8. https://CRAN.R-project.org/package=tibble
##   - Kuznetsova A, Brockhoff PB, Christensen RHB (2017). "lmerTest Package:Tests in Linear Mixed Effects Models." _Journal of StatisticalSoftware_, *82*(13), 1-26. doi: 10.18637/jss.v082.i13 (URL:https://doi.org/10.18637/jss.v082.i13).
##   - Lionel Henry and Hadley Wickham (2020). purrr: Functional Programming Tools. R package version 0.3.4. https://CRAN.R-project.org/package=purrr
##   - Lüdecke D (2018). "ggeffects: Tidy Data Frames of Marginal Effects fromRegression Models." _Journal of Open Source Software_, *3*(26), 772.doi: 10.21105/joss.00772 (URL: https://doi.org/10.21105/joss.00772).
##   - R Core Team (2021). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.
##   - Rinker, T. W. & Kurkiewicz, D. (2017). pacman: Package Management for R. version 0.5.0. Buffalo, New York. http://github.com/trinker/pacman
##   - Wickham et al., (2019). Welcome to the tidyverse. Journal of Open Source Software, 4(43), 1686, https://doi.org/10.21105/joss.01686
##   - Yihui Xie (2021). knitr: A General-Purpose Package for Dynamic Report Generation in R. R package version 1.37.